Numpy 介紹

In [1]:
# 起手式
import numpy as np

建立 ndarray

array 裡面必須是同一種資料型態

In [32]:
a = np.array([1,2,3,4]); print(a, repr(a))
b = np.array([1,2,3,'abc']); print(b, repr(b))
c = np.array([1,2,3,'abc', ['cde', 123]]); print(c)

[1 2 3 4] array([1, 2, 3, 4])
['1' '2' '3' 'abc'] array(['1', '2', '3', 'abc'], 
      dtype='<U21')

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-32-583021b3090d> in <module>()
      1 a = np.array([1,2,3,4]); print(a, repr(a))
      2 b = np.array([1,2,3,'abc']); print(b, repr(b))
----> 3 c = np.array([1,2,3,'abc', ['cde', 123]]); print(c)

ValueError: setting an array element with a sequence
In [35]:
a[2]
Out[35]:
3
In [10]:
a[2] = 999; a[2]
Out[10]:
999
In [44]:
x = np.array([1,2,3,4])
In [45]:
y = np.array([[1.,2,3],[4,5,6]])
y
Out[45]:
array([[ 1.,  2.,  3.],
       [ 4.,  5.,  6.]])
In [46]:
y[1, 0]
Out[46]:
4.0

看 ndarray 的第一件事情： shape , dtype

In [47]:
x.shape
Out[47]:
(4,)
In [48]:
y.shape # 原始建構 list 的層級由外往內的 dimension
Out[48]:
(2, 3)
In [49]:
x.dtype
Out[49]:
dtype('int64')
In [50]:
y.dtype
Out[50]:
dtype('float64')

有時候，可以看圖

In [132]:
# import matplotlib
%matplotlib inline
from matplotlib import pyplot as plt
plt.style.use('ggplot')

# 畫圖
plt.plot(x, 'x');

有很多其他建立的方式

In [10]:
# 建立 0 array
np.zeros_like(y)
Out[10]:
array([[ 0.,  0.,  0.],
       [ 0.,  0.,  0.]])
In [40]:
np.zeros((10,10))
Out[40]:
array([[ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.]])
In [43]:
np.ones((3, 3, 2))
Out[43]:
array([[[ 1.,  1.],
        [ 1.,  1.],
        [ 1.,  1.]],

       [[ 1.,  1.],
        [ 1.,  1.],
        [ 1.,  1.]],

       [[ 1.,  1.],
        [ 1.,  1.],
        [ 1.,  1.]]])
In [98]:
# 跟 range 差不多
x = np.arange(0, 10, 0.1); x
Out[98]:
array([ 0. ,  0.1,  0.2,  0.3,  0.4,  0.5,  0.6,  0.7,  0.8,  0.9,  1. ,
        1.1,  1.2,  1.3,  1.4,  1.5,  1.6,  1.7,  1.8,  1.9,  2. ,  2.1,
        2.2,  2.3,  2.4,  2.5,  2.6,  2.7,  2.8,  2.9,  3. ,  3.1,  3.2,
        3.3,  3.4,  3.5,  3.6,  3.7,  3.8,  3.9,  4. ,  4.1,  4.2,  4.3,
        4.4,  4.5,  4.6,  4.7,  4.8,  4.9,  5. ,  5.1,  5.2,  5.3,  5.4,
        5.5,  5.6,  5.7,  5.8,  5.9,  6. ,  6.1,  6.2,  6.3,  6.4,  6.5,
        6.6,  6.7,  6.8,  6.9,  7. ,  7.1,  7.2,  7.3,  7.4,  7.5,  7.6,
        7.7,  7.8,  7.9,  8. ,  8.1,  8.2,  8.3,  8.4,  8.5,  8.6,  8.7,
        8.8,  8.9,  9. ,  9.1,  9.2,  9.3,  9.4,  9.5,  9.6,  9.7,  9.8,
        9.9])
In [99]:
# 亂數
y = np.random.uniform(-1,1, size=x.shape)
plt.plot(x, y)
Out[99]:
[<matplotlib.lines.Line2D at 0x1088cddd8>]
In [117]:
# 亂數
x = np.arange(0, 10, 0.1)
y = np.random.normal(-1,1, size=x.shape); print(x.shape)
print(y)

plt.figure(figsize=(10,4))
plt.subplot(121)
plt.plot(x, y)

plt.subplot(122)
plt.plot(np.ones_like(y), y, 'o')

(100,)
[ -1.77536057e+00  -5.58850085e-01   1.23603129e-01  -2.79965516e+00
  -1.16043732e-01  -3.31023913e+00   1.40106984e-01   1.25954644e+00
  -1.61128993e-02  -2.22375287e+00  -1.52121791e+00  -1.25577096e+00
   1.01228894e+00   3.27826233e-01  -1.48552766e+00  -1.80007049e+00
  -2.07331247e+00   5.19667031e-02  -8.91206573e-01   8.41519462e-01
  -1.86567280e+00  -4.30903941e+00  -1.94357447e+00  -1.50376068e+00
  -1.38165144e+00  -1.77447778e+00  -1.71507809e+00  -8.23994137e-01
  -1.50939456e-01  -2.86144661e+00  -1.30062976e+00  -1.16461036e+00
   7.96455741e-01  -9.59019031e-01  -8.12294330e-01   1.07209799e-02
  -1.10333505e+00  -1.81039207e+00   9.59588702e-02  -1.95571024e-01
   5.77956065e-01  -1.50995679e+00  -1.70879121e-01  -1.77612337e+00
  -4.58137850e-01  -6.33671594e-01  -1.01754735e+00  -8.05856148e-01
  -1.92429979e+00  -2.16236466e+00  -3.29536520e-01  -1.19372061e+00
  -2.34612701e+00  -4.19097869e-01  -4.67284003e-01  -2.10393553e-01
   3.47383667e-01  -1.02380581e+00  -3.77506211e-01  -1.92873295e+00
   5.40820417e-01  -7.12943121e-01  -1.36345358e+00  -1.71962459e+00
  -2.38928257e-01  -8.86216279e-01  -8.02588461e-01  -2.89091353e+00
  -1.61277592e+00  -7.62655039e-01   3.61217144e-01  -1.14167118e+00
  -7.98138623e-01  -8.81683323e-01  -1.92054753e-01  -9.23259200e-04
  -2.80090776e+00  -4.72684068e-01   4.63087714e-02  -1.17330090e+00
  -2.13574058e+00  -6.98669027e-01   3.19982674e-01  -1.85260007e+00
  -1.09808729e+00  -7.52176246e-01  -1.71872906e+00  -2.24650350e-01
  -1.51593860e+00  -8.35140790e-01  -1.91006103e+00  -9.20731598e-01
  -1.22140403e+00   6.21012317e-01  -7.56384178e-01  -1.83435263e+00
  -8.63687547e-01  -1.08539381e+00  -1.41929657e+00  -1.22004071e+00]
Out[117]:
[<matplotlib.lines.Line2D at 0x108fb8828>]

這是一堆資料

  • 資料有什麼資訊？
  • 資料有什麼限制？
  • 這些限制有什麼意義？好處？
  • 以前碰過什麼類似的東西？
  • 可以套用在哪些東西上面？
  • 可以怎麼用（運算）?

最簡單的計算是 逐項計算

see also np.vectorize

In [350]:
np.array([1, 2, 3]) * np.array([100, 200])

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
/Users/leoluyi/Dropbox/Courses/PyCon2017_Numpy_tutorial/numpy_tutorial/numpy/q_scale2.py in <module>()
----> 1 np.array([1, 2, 3]) * np.array([100, 200])

ValueError: operands could not be broadcast together with shapes (3,) (2,)
In [351]:
np.array([1, 2, 3]) * np.array([100, 200, 300])
Out[351]:
array([100, 400, 900])
In [101]:
x = np.linspace(0, 2* np.pi, 1000) # 0~2pi 平均切 1000 個點
plt.plot(x, np.sin(x))
Out[101]:
[<matplotlib.lines.Line2D at 0x1083c62e8>]

Q0:

畫出 $y=x^2+1$ 或其他函數的圖形 用

plt.plot?

看看 plot 還有什麼參數可以玩

In [130]:
x = np.linspace(-2, 2, 30)
y = x**2 + 1
print(plt.style.available)
with plt.style.context((['fivethirtyeight'])):
    plt.plot(x, y, '', x, y, 'o')

['classic', 'fivethirtyeight', 'seaborn-talk', 'dark_background', 'seaborn-muted', 'seaborn-pastel', 'ggplot', 'bmh', 'seaborn-notebook', 'seaborn-bright', 'seaborn-poster', 'seaborn-paper', 'seaborn-deep', 'seaborn-dark-palette', 'seaborn-colorblind', 'seaborn-dark', 'seaborn-ticks', 'seaborn-whitegrid', 'seaborn-darkgrid', 'seaborn-white', 'grayscale']
In [134]:
# 或者看看參考範例
# %run q0.py
# %load q0.py
x = np.linspace(-1,1,11)
plt.plot(x, x**2+1, "g^", x, x**3, 'ro', x, x**2+1, 'b');

Q1:

試試看圖片。 使用

from PIL import Image
# 讀入 PIL Image (這張圖是從 openclipart 來的 cc0)
img = Image.open('img/Green-Rolling-Hills-Landscape-800px.png')
# 圖片轉成 ndarray
img_array = np.array(img)
# ndarray 轉成 PIL Image
Image.fromarray(img_array)

看看這個圖片的內容, dtype 和 shape

In [139]:
# 參考答案
# %load q1.py
from PIL import Image
img = Image.open('img/Green-Rolling-Hills-Landscape-800px.png')
img_array = np.array(img)
print("img_array.shape={}".format(img_array.shape))
print("img_array.dtype={}".format(img_array.dtype))
Image.fromarray(img_array)

img_array.shape=(530, 800, 4)
img_array.dtype=uint8
Out[139]:

Indexing

可以用類似 list 的 indexing

In [17]:
a = np.arange(30)
a
Out[17]:
array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12, 13, 14, 15, 16,
       17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29])
In [18]:
a[5]
Out[18]:
5
In [19]:
a[3:7]
Out[19]:
array([3, 4, 5, 6])
In [20]:
# 列出所有奇數項
a[1::2]
Out[20]:
array([ 1,  3,  5,  7,  9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29])
In [21]:
# 還可以用來設定值
a[1::2]  = -1
a
Out[21]:
array([ 0, -1,  2, -1,  4, -1,  6, -1,  8, -1, 10, -1, 12, -1, 14, -1, 16,
       -1, 18, -1, 20, -1, 22, -1, 24, -1, 26, -1, 28, -1])
In [22]:
# 或是
a[1::2] = -a[::2]-1
a
Out[22]:
array([  0,  -1,   2,  -3,   4,  -5,   6,  -7,   8,  -9,  10, -11,  12,
       -13,  14, -15,  16, -17,  18, -19,  20, -21,  22, -23,  24, -25,
        26, -27,  28, -29])

Q2

給定

x = np.arange(30)
a = np.arange(30)
a[1::2] = -a[1::2]

畫出下面的圖

In [140]:
%run -i q2.py
#%load q2.py
In [171]:
x = np.arange(30); print(x.shape)
a = np.arange(30)
a[1::2] = -a[1::2]
px = x[0::2]
py = x[0::2]; print(py.shape)
qx = x[1::2]
qy = -x[1::2]; print(qy.shape)

with plt.style.context('bmh'):
    plt.plot(x, a, '', px, py, 'x', qx, qy, 'v')

(30,)
(15,)
(15,)

ndarray 也可以

In [199]:
b = np.array([[1,2,3], [4,5,6], [7,8,9]])
b
Out[199]:
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
In [200]:
b[1][2]
Out[200]:
6
In [201]:
b[1,2] # 1,2 本身就是 tuple
Out[201]:
6
In [180]:
b[(1,2)]
Out[180]:
6
In [179]:
b[1]
Out[179]:
array([4, 5, 6])
In [187]:
b[1,:]
Out[187]:
array([4, 5, 6])
In [196]:
from pprint import pprint
pprint(img_array[:,:, 2]) # blue channel
Image.fromarray(img_array[:,:, 2])

array([[242, 241, 232, ..., 238, 239, 239],
       [243, 243, 235, ..., 239, 239, 239],
       [243, 243, 238, ..., 238, 239, 240],
       ..., 
       [  3,   2,   3, ...,   4,   4,   3],
       [  3,   0,  34, ...,   4,   4,   4],
       [  3,   0,   7, ...,   4,   7,  15]], dtype=uint8)
Out[196]:

Q3

動手試試看各種情況 比方

b = np.random.randint(0,99, size=(10,10))
b[::2, 2]
In [220]:
b = np.random.randint(0,99, size=(10,10))
print(b[::2, 2])
# slice 物件
print(b[slice(0, 10, 2), 2])

[84 47 33 98 43]
[84 47 33 98 43]

Fancy indexing

物件會 copying

In [344]:
np.random.seed(1)
b = np.random.randint(0,99, size=(5,10))
b
Out[344]:
array([[37, 12, 72,  9, 75,  5, 79, 64, 16,  1],
       [76, 71,  6, 25, 50, 20, 18, 84, 11, 28],
       [29, 14, 50, 68, 87, 87, 94, 96, 86, 13],
       [ 9,  7, 63, 61, 22, 57,  1,  0, 60, 81],
       [ 8, 88, 13, 47, 72, 30, 71,  3, 70, 21]])

試試看下面的結果

想一下是怎麼一回事(numpy 在想什麼？)

In [334]:
b[[1,3]] #  == b[([1,3],)] 利用 list 做 indexing
Out[334]:
array([[76, 71,  6, 25, 50, 20, 18, 84, 11, 28],
       [ 9,  7, 63, 61, 22, 57,  1,  0, 60, 81]])
In [335]:
np.array_equal( b[[1,3]], b[([1,3],)] )
Out[335]:
True
In [336]:
b[(1,3)] # == b[1,3]
Out[336]:
25
In [337]:
b[[1,2], 
  [3,4]]  # 先篩 d1 的 1,2，再篩 d2 的 3,4
Out[337]:
array([25, 87])
In [338]:
b[[(1,2),
   (3,4)]] # 同上，先篩 d1 的 1,2，再篩 d2 的 3,4
Out[338]:
array([25, 87])
In [339]:
b[[1,2]][:,[3,4]] # 同時篩 d1: (1, 3), d2: (3,4)
Out[339]:
array([[25, 50],
       [68, 87]])
In [340]:
b[[True, False, False, True, False]]
Out[340]:
array([[37, 12, 72,  9, 75,  5, 79, 64, 16,  1],
       [ 9,  7, 63, 61, 22, 57,  1,  0, 60, 81]])
In [345]:
print(b)
b[:, np.arange(10) % 3 == 0]

[[37 12 72  9 75  5 79 64 16  1]
 [76 71  6 25 50 20 18 84 11 28]
 [29 14 50 68 87 87 94 96 86 13]
 [ 9  7 63 61 22 57  1  0 60 81]
 [ 8 88 13 47 72 30 71  3 70 21]]
Out[345]:
array([[37,  9, 79,  1],
       [76, 25, 18, 28],
       [29, 68, 94, 13],
       [ 9, 61,  1, 81],
       [ 8, 47, 71, 21]])

Q4

把 b 中的偶數都變成 -1

In [34]:
#參考範例
%run -i q4.py

array([[13, -1, 41, 43, -1, -1, 93, 55, -1, -1],
       [89, -1, 15, 21, -1, -1, 65, -1, -1,  7],
       [-1, -1, -1, 25, -1, -1, -1, 31, -1, -1],
       [-1, -1, -1, 97, 45, -1,  9, -1, 15, -1],
       [73, -1, -1, -1, 33, 45, -1, -1, -1, -1]])
In [328]:
b = np.random.randint(0,99, size=(5,10)); print('id of b: {}'.format(id(b)))
pprint(b)
b[b % 2 == 0] = -1; print('id of b: {}'.format(id(b)))
pprint(b)

id of b: 4506298768
array([[86, 70, 66, 71, 48, 54, 15,  5, 17, 42],
       [20, 48, 22, 13, 97, 53, 84, 10, 96, 55],
       [61, 56, 89, 21, 96, 83, 25, 14, 13, 84],
       [43,  6, 77, 56, 59, 15, 24,  9, 66, 71],
       [53, 69, 36, 21, 40, 77, 91, 49, 47, 77]])
id of b: 4506298768
array([[-1, -1, -1, 71, -1, -1, 15,  5, 17, -1],
       [-1, -1, -1, 13, 97, 53, -1, -1, -1, 55],
       [61, -1, 89, 21, -1, 83, 25, -1, 13, -1],
       [43, -1, 77, -1, 59, 15, -1,  9, -1, 71],
       [53, 69, -1, 21, -1, 77, 91, 49, 47, 77]])

用圖形來練習

In [35]:
# 還記得剛才的
from PIL import Image
img = Image.open('img/Green-Rolling-Hills-Landscape-800px.png')
img_array = np.array(img)
Image.fromarray(img_array)
Out[35]:
In [235]:
# 用來顯示圖片的函數
from IPython.display import display
def show(img_array):
    display(Image.fromarray(img_array))

Q. Scaling

  • 將圖片縮小成一半
  • 擷取中間一小塊
  • 圖片上下顛倒
  • 左右鏡射
  • 去掉綠色
  • 將圖片放大兩倍
  • 貼另外一張圖到大圖中
    from urllib.request import urlopen
url = "https://raw.githubusercontent.com/playcanvas/engine/master/examples/images/animation.png"
simg = Image.open(urlopen(url))
  • 紅綠交換
  • 團片變成黑白 參考 Y=0.299R+0.587G+0.114B
    • 會碰到什麼困難？ 要如何解決
In [37]:
# 將圖片縮小成一半
%run -i q_half.py
In [274]:
# 將圖片放大
# %run -i q_scale2.py
from urllib.request import urlopen
url = "https://raw.githubusercontent.com/playcanvas/engine/master/examples/images/animation.png"
simg = Image.open(urlopen(url))
simg_array = np.array(simg)
print("原始大小")
show(simg_array)

print("增高兩倍")
simg_array_h2 = simg_array[[i//2 for i in range(2*simg_array.shape[0])],:]
show(simg_array_h2)

print("增寬兩倍")
simg_array_w2 = simg_array[:, [i//2 for i in range(2*simg_array.shape[1])]]
show(simg_array_w2)

print("放大兩倍")
x_idx = [i//2 for i in range(2*simg_array.shape[0])]
y_idx = [i//2 for i in range(2*simg_array.shape[1])]
simg_array_2 = simg_array[x_idx][:,y_idx]

show(simg_array_2)

原始大小

增高兩倍

增寬兩倍

放大兩倍
In [236]:
# 圖片上下顛倒
show(img_array[::-1])
In [365]:
# %run -i q_paste.py
simg_array = np.array(simg); #print(simg_array)
img_array2 = img_array.copy()
print("簡單的")
img_array2[200:400, 300:500] = simg_array
show(img_array2)

print("這樣呢？")
img_array2 = img_array.copy()
# print(np.apply_along_axis(lambda a: np.sum(a[:3]), 2, simg_array))
is_black = np.apply_along_axis(lambda a: np.sum(a[:3]), 2, simg_array) > 10; print(is_black)
img_array2[200:400, 300:500][is_black] = simg_array[is_black]
show(img_array2)

簡單的

這樣呢？
[[False False False ..., False False False]
 [False False False ..., False False False]
 [False False False ..., False False False]
 ..., 
 [False False False ..., False False False]
 [False False False ..., False False False]
 [False False False ..., False False False]]
In [41]:
%run -i q_grayscale.py

Q

  • 挖掉個圓圈？ (300,300)中心，半徑 100
  • 旋轉九十度? x,y 互換？
In [42]:
# 用迴圈畫圓
%run -i q_slow_circle.py
In [43]:
# 用 fancy index 畫圓
%run -i q_fast_circle.py
In [395]:
x = np.arange(800).reshape(1, -1); print(x.shape)
y = np.arange(530).reshape(-1, 1); print(y.shape)
print((x + y).shape) # broadcasting
in_circle = (x-400)**2 + (y-100)**2 <- 80**2 # broadcasting circle

(1, 800)
(530, 1)
(530, 800)
Out[395]:
array([[False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False],
       ..., 
       [False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False],
       [False, False, False, ..., False, False, False]], dtype=bool)

indexing 的其他用法

In [44]:
# 還可以做模糊化
a = img_array.astype(float)
for i in range(10):
    a[1:,1:] = (a[1:,1:]+a[:-1,1:]+a[1:,:-1]+a[:-1,:-1])/4
show(a.astype('uint8'))
In [45]:
# 求邊界
a = img_array.astype(float)
a = a @ [0.299, 0.587, 0.114, 0]
a = np.abs((a[1:]-a[:-1]))*2
show(a.astype('uint8'))

Reshaping

.flatten 拉平看看資料在電腦中如何儲存？

查看 .reshape, .T, np.rot00, .swapaxes .rollaxis 然後再做一下上面的事情

  • .T： 最後兩軸交換
  • swapaxes: 交換軸
  • moveaxis: 移動軸位置
  • rollaxis: 軸往前往後移動
In [377]:
A = np.random.randint(0, 10, size=(6,6))
A
Out[377]:
array([[1, 2, 7, 2, 6, 0],
       [9, 2, 6, 6, 2, 7],
       [7, 0, 6, 5, 1, 4],
       [6, 0, 6, 5, 1, 2],
       [1, 5, 4, 0, 7, 8],
       [9, 5, 7, 0, 9, 3]])
In [381]:
print(A.flatten()) # copy
A.ravel()
print(A)

[1 2 7 2 6 0 9 2 6 6 2 7 7 0 6 5 1 4 6 0 6 5 1 2 1 5 4 0 7 8 9 5 7 0 9 3]
[[1 2 7 2 6 0]
 [9 2 6 6 2 7]
 [7 0 6 5 1 4]
 [6 0 6 5 1 2]
 [1 5 4 0 7 8]
 [9 5 7 0 9 3]]
In [382]:
A.reshape((2,2,3,-1))
Out[382]:
array([[[[1, 2, 7],
         [2, 6, 0],
         [9, 2, 6]],

        [[6, 2, 7],
         [7, 0, 6],
         [5, 1, 4]]],


       [[[6, 0, 6],
         [5, 1, 2],
         [1, 5, 4]],

        [[0, 7, 8],
         [9, 5, 7],
         [0, 9, 3]]]])
In [400]:
# reshaping 的應用
R,G,B,A = img_array.reshape(-1,4).T
plt.hist((R,G,B,A), color="rgbk");

堆疊在一起

查看:

  • np.vstack
  • np.hstack
  • np.concatenate

然後試試看

In [47]:
# 例子
show(np.hstack([img_array, img_array2]))
In [48]:
# 例子
np.concatenate([img_array, img_array2], axis=2).shape
Out[48]:
(530, 800, 8)

作用在整個 array/axis 的函數

In [49]:
np.max([1,2,3,4])
Out[49]:
4
In [50]:
np.sum([1,2,3,4])
Out[50]:
10
In [51]:
np.mean([1,2,3,4])
Out[51]:
2.5
In [52]:
np.min([1,2,3,4])
Out[52]:
1
In [402]:
np.argmax([2,1,5,4])
Out[402]:
2

多重意義的運用， 水平平均，整合垂直平均

In [53]:
x_mean = img_array.astype(float).min(axis=0, keepdims=True)
print(x_mean.dtype, x_mean.shape)
y_mean = img_array.astype(float).min(axis=1, keepdims=True)
print(y_mean.dtype, y_mean.shape)
# 自動 broadcast 
xy_combined = ((x_mean+y_mean)/2).astype('uint8')
show(xy_combined)

float64 (1, 800, 4)
float64 (530, 1, 4)

Tensor 乘法

In [414]:
x = np.zeros((3,4,5))
print(x[None, :, :, :].shape) # 增維
print(x[:, None, :, :].shape) # 增維
print(x[:, None, ..., np.newaxis].shape) # 增維

(1, 3, 4, 5)
(3, 1, 4, 5)
(3, 1, 4, 5, 1)

先從點積開始

In [54]:
# = 1*4 + 2*5 + 4*6
np.dot([1,2,3], [4,5,6])
Out[54]:
32
In [428]:
u=np.array([1,2,3])
v=np.array([4,5,6])
print( u@v ) # 矩陣乘法 np.matmult
print(np.matmul(u,v))
print( (u*v).sum() )

32
32
32

矩陣乘法

如果忘記矩陣乘法是什麼了， 參考這裡 http://matrixmultiplication.xyz/ 或者 http://eli.thegreenplace.net/2015/visualizing-matrix-multiplication-as-a-linear-combination/

矩陣乘法可以看成是:

  • 所有組合（其他軸）的內積（共有軸）
  • 多個行向量線性組合
  • 代入線性方程式 A1-矩陣與基本列運算.ipynb
  • 用 numpy 來理解
    np.sum(a[:,:, np.newaxis] * b[np.newaxis, : , :], axis=1)
dot(a, b)[i,k] = sum(a[i,:] * b[:, k])

高維度

要如何推廣？

  • tensordot, tensor contraction, a.shape=(3,4,5), b.shape=(4,5,6), axis = 2 時等價於
np.sum(a[..., np.newaxis] * b[np.newaxis, ...], axis=(1, 2))
tensordot(a,b)[i,k]=sum(a[i, ...]* b[..., k])

https://en.wikipedia.org/wiki/Tensor_contraction

  • dot
dot(a, b)[i,j,k,m] = sum(a[i,j,:] * b[k,:,m])
  • tensordot
np.tensordot(a,b, axes=(-1,-2))
  • matmul 最後兩個 index 當成 matrix
a=np.random.random(size=(3,4,5))
b=np.random.random(size=(3,5,7))
(a @ b).shape
np.sum(a[..., np.newaxis] * np.moveaxis(b[..., np.newaxis], -1,-3), axis=-2)

np.einsum('ii', a) # trace(a)
np.einsum('ii->i', a) #diag(a)
np.einsum('ijk,jkl', a, b) # tensordot(a,b)
np.einsum('ijk,ikl->ijl', a,b ) # matmul(a,b)
In [416]:
A=np.random.randint(0,10, size=(5,3))
A
Out[416]:
array([[9, 1, 4],
       [4, 6, 8],
       [8, 9, 2],
       [7, 5, 5],
       [4, 5, 8]])
In [417]:
B=np.random.randint(0,10, size=(3,7))
B
Out[417]:
array([[5, 8, 1, 1, 8, 7, 0],
       [3, 4, 2, 0, 3, 5, 1],
       [2, 4, 3, 0, 6, 0, 7]])
In [419]:
A.dot(B).shape
Out[419]:
(5, 7)

Q

  • 手動算算看 A,B 的 dot
  • 試試看其他的乘法
In [418]:
np.matmul(A, B).shape
Out[418]:
(5, 7)
In [423]:
np.multiply(A, A)
Out[423]:
array([[81,  1, 16],
       [16, 36, 64],
       [64, 81,  4],
       [49, 25, 25],
       [16, 25, 64]])
In [447]:
C = np.random.randint(0, 10, size=(3,5,6))
D = np.random.randint(0, 10, size=(5,6,7))
np.tensordot(C, D, axes=2).shape # (default) tensor double contraction
Out[447]:
(3, 7)

小結

numpy 以 ndarray 為中心

  • 最基本的運算是逐項運算(用 np.vectorize把一般的函數變成逐項運算 )
  • indexing 很好用
  • reshaping
  • 整合的操作與計算